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This paper describes the development of a numerical model that couples the thermal 
interaction between ice particles, water droplets, and the flowing gas of an icing wind tunnel 
for simulation of NASA Glenn Research Center’s Propulsion Systems Laboratory (PSL). The 
ultimate goal of the model is to better understand the complex interactions between the test 
parameters and have greater confidence in the conditions at the test section of the PSL tunnel. 
The model attempts to explain the observed changes in test conditions by coupling the 
conservation of mass and energy equations for both the cloud particles and flowing gas mass. 
The model uses isentropic relations to relate gas temperature, velocity, density and pressure 
with respect to the PSL geometry. Measurements were taken at the PSL during wind tunnel 
tests simulating ice-crystal and mixed-phase icing that relate to ice accretions within turbofan 
engines in May 2015. The model was compared to experimentally measured values, where test 
conditions varied gas temperature, pressure, velocity and humidity levels, as well as the cloud 
total water content, particle initial temperature, and particle size distribution. Wet-bulb 
temperatures were generally within a few degrees of freezing. The model showed good 
agreement with experimentally measured values, to within approximately 30% of the 
measured change in gas temperature and humidity at the tunnel test section. The model did 
reasonably well in predicting melt content (liquid mass to total mass) at the test section, 
especially for clouds with larger particle sizes. In addition, the model predicted particle size 
at the tunnel exit with good agreement, however, the comparison was limited to clouds 
consisting of a small particle size distribution. One of the key findings from this work is that 
there was a nearly constant but slight increase in total wet-bulb temperature when the spray 
cloud was activated for every test and simulation. In addition, the total wet-bulb temperature 
in the tunnel plenum was a large factor in determining cloud phase. 


Nomenclature 

area (m7) 
= speed of sound (m s"!) 
= specific heat capacity (J kg! K") 
Cloud Droplet Probe (by Droplet Measurement Technology, Inc.) 
= coefficient of drag 
= particle diameter (m or pm) 
= volumetric diameter (m or um) 
= enthalpy rate (W) 

= heat transfer coefficient (W m? K"!) 
1 = mass transfer coefficient (m s*!) 
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Subscripts 
air = 


bulk = 


Yodiff = 


NASA Icing Research Tunnel 

ice water content (g m*) 

ratio of specific heats of a gas 

kinetic energy rate (W) 

latent heat of phase change (J kg"!) 

liquid water content (g m=) 

median volume diameter (m or 1m) 

mass (kg) 

National Research Council of Canada 
pressure (Pa) 

particle size distribution 

NASA Propulsion Systems Laboratory 
heat flux (W) 

NRC’s Research Altitude Test Facility 
relative humidity (%) 

temperature or gas temperature (K, °C, or °F) 
test point # 

wet-bulb temperature (K, °C, or °F) 

total water content (g m°) 

time (s) 

velocity (m s‘!) 

distance (m) 

change 

differential 

melt fraction (liquid water mass/total water mass) 
density (kg m’) 

mixed mass ratio (vapor mass/dry air mass) 
total number of particles in a bin 


dry air 

reference to bulk total water content 
convection 

exit 

experiment 

gas 

inlet condition or i number 
in (into a control volume) 
melt 

number of bins 

nozzle 

spray on 

spray off 

out (out of a control volume) 
particle (ice, water, or mixed phase) 
reference 

static 

simulation 

water 

water vapor 

total or plenum 

percent difference 
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I. Introduction 


EVERAL jet engine power-loss events of commercial aircraft at high altitude have been reported since the 1990’s. 

Mason et al.! discussed how power-loss could result from ice crystals entering the engine core, partially melting 
from the warm engine gas flow and refreezing on internal components. To better understand this phenomenon and 
determine the physical mechanism of icing in the core of the engine, NASA performed fundamental physics of 
ice-crystal ice accretion tests. A two-week test effort was conducted at the NASA Propulsion Systems Laboratory 
(PSL) in March 2016 with preliminary data collected in May 2015. During these tests it was critical to quantify key 
icing parameters at the test section, such as temperature, humidity, pressure, velocity, total water content, fraction of 
liquid to total water content, and the particle size distribution. The main objective of these tests was to generate a 
prescribed mixed-phase cloud at the test section. 

During the PSL fundamental physics of ice-crystal icing tests, and on previous occasions”™, it has been observed 
that the test conditions, most notably temperature and humidity, change when the icing cloud is activated. In some 
cases, the gas temperature decrease was considerable, by as much as several degrees Celsius, accompanied by a 
measurable increase in water vapor. It is hypothesized that the ice particles and water droplets thermally interact with 
the flowing gas causing the gas temperature and humidity to change by the time the cloud particles reache the test 
section. 

Existing thermal models®"' show little effect due to coupling of the ice/water particles with the flowing gas, or do 
not couple, approximating the gas mass as an infinite thermal reservoir with unchanging properties. A thermal model, 
previously written by Bartkus et al.!? attempts to explain the observed changes in test conditions by coupling the 
conservation of mass and energy equations for both the ice/water particle and flowing gas mass. The model simulates 
both ice particle melting and droplet freezing, including spontaneous latent heat release during freezing for 
supercooled droplets. In the model, humidity changes are due to ice particle sublimation/deposition and water droplet 
evaporation/condensation. The associated latent heating/cooling and mass transfer change the ice water content (WC) 
and liquid water content (LWC) from the moment of injection to the time the flowing mass reaches the test section. 
As a result, the model predicts changes in gas temperature, humidity, WC and LWC, as well as the particle phase, 
size, temperature and velocity at the test section. To more accurately simulate spray injection, the model incorporates 
a distribution of particle sizes (24+ bins) as an initial input condition. 

The model was specifically written to compare experimental data measured during a collaboratively conducted 
icing campaign between NASA GRC and the National Research Council of Canada (NRC) at the Research Altitude 
Test Facility (RATFac)*°. The model generally predicted the correct trend, but did not fully explain the changes that 
were observed at the RATFac icing wind tunnel. In general, the model was able to account for approximately 20% of 
the water that vaporized from the icing cloud, and approximately 20% of the gas temperature change that was observed 
experimentally!” One of the reasons offered for this discrepancy may be that water/ice film on the tunnel walls that 
had been observed during testing — which can increase vapor content and temperature change — was not accounted for 
in the model. 

The NASA PSL is an altitude jet-engine test facility that generates ice particles using a liquid water spray nozzle 
system whereby the injected water droplets can freeze out as the cloud flows towards the engine/test section’: 4. The 
droplets freeze due to a combination of convective heat transfer and evaporative cooling. Unlike at RATFac, PSL can 
chose to not to deposit ice or water on the tunnel walls by turning off individual spray nozzles. This allows for a better 
platform for comparing the present thermal model — with modifications specific to the PSL — to experimentally 
measured data. The model has currently been modified to include isentropic relations with respect to the PSL 
geometry, which is described later in the paper. With this modification, the model aided in finding a range of conditions 
that resulted in a mixed-phase cloud at the exit of the PSL free jet duct during the 2015 and 2016 tests. Data measured 
during the recent March 2016 tests must still be analyzed, therefore this paper will only compare the data collected 
during preliminary testing in May 2015 with predicted model results. A detailed description of the preliminary 2015 
testing, and an explanation of the data collected can be found elsewhere!> '°, 
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II. Thermal Model Description 


The thermal model simulates an icing wind tunnel by calculating particle and gas properties using expressions 
derived from the conservation of mass, momentum and energy equations. The particle mass and gas mass equations 
were fully coupled, as were the particle energy and gas energy equations. The conservation of momentum equation is 
solved only in reference to the particle, as it is approximated that the particles’ effect on the flowing gas is negligible. 
Differential equations for change in gas temperature (Teas), particle temperature (7,), gas mass (iMgas), particle mass 
(mp), gas velocity (Vegas), particle velocity (v,), particle melt ratio (7), gas density (Peas), and pressure (P) were derived 
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from equations of conservation. Since the tunnel geometry changes axially, isentropic relations were incorporated into 
the model to relate gas temperature, velocity, density and pressure with respect to changing are of the PSL geometry. 


A. Assumptions 
The thermal model operates under the following assumptions: 


B. Ther 


The values along the tunnel axis are steady (not a transient problem). 

Gas and particle flow are steady and one dimensional. 

The dry air and water vapor (that composes the gas) are treated as ideal gases. 

Gas is well mixed, meaning that temperature, water vapor content and all thermodynamic properties are 
homogeneous within the gas control volume at any axial location (or time). 

The system (tunnel) is adiabatic and mass is conserved. 

Particle size distribution is characterized by a discrete set of diameters. 

Particles of every size in the distribution are evenly spaced within the gas from the injection point to the test 
section (i.e. uniform total water content for a given cross-section). 

All particles are perfectly spherical. 

Particle aggregation/coalescences and breakup through collision are negligible. 

Particles are injected in the direction of the flow and remain entrained with negligible gravity affects. 
Temperature is uniform within the particle. 

Mixed phase particles are spatially homogeneous in water/ice content. This means that evaporation and 
sublimation occurs from a mixed phase particle surface and the rate is determined by the water/ice content 
ratio. This applies to the reverse process of condensation and deposition as well. 

Evaporation, sublimation, condensation and deposition occur at the particle surface at particle temperature. 
The flow of particles and gas is a continuous stream. This means that while following a particular set of 
particles within a reference control volume, faster moving particles and gas that exit the reference control 
volume are replaced by thermodynamically identical particles and gas from a neighboring control volume 
that is upstream of the reference control volume. The reference control volume is centered on the slowest 
moving particle in a distribution of particles. 


mal Model Equations — Gas Conservation Equations 


The changing geometry of the tunnel must be addressed in the conservation of gas equations. The equations to 
follow refer to the infintesimal control volume in Figure 1. The volume is a representative of an infinetly thin cross- 


section o 


f the tunnel. The flow is approximated to be one dimensional, flowing from left to right in Figure 1. With 


respect to nomenclature in this section, a parameter with a subscript of ‘gas’, refers to a property of the flowing gas 
(combined air and vapor), while the subscript, ‘wv’, refers to a property of water vapor only. Also, the ‘dot’ accent, 
for example m, refers to a rate. 


Flow 


Direction 
eh 


Figure 1: Differential control volume of width ox. 


A unit mass of gas, mgas, enters the control volume 
with a velocity of Vegas, a temperature of Tzas, a density 
— ox _ of Pgas, and a pressure of P. Water vapor mass, Mwy, 

enters into the control volume as a mass source (or a 
sink if vapor is condensing), with velocity vy, and 
temperature of 7,,,. The gas and vapor masses enter 
the contol volume with an inlet surface area of A. At 
an infintesimal distance downstream, Ox, the vapor 
mass exits the volume with velocity Vgas OVeas, a 
temperature of Tgast+OT gas, a density of PeastOPeas, and 
a pressure of P+0P. The gas exits the control volume 
with an exit surface area of A+0A. 

Expressions that describe their changes through 
the infinitesimal volume for the five mentioned 
variables (Vegas, Teas, Pgas, P, and A) are desired. Four 
of the five variables can be expressed using the 
conservation of mass, momentum and _ energy 
equations, along with the speed of sound equation. 
The fifth variable, the change in area, can be 
expressed using the known geomety. 


Veas 
Tyas 
Peas 
P 

A 
Meas 
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The mass flow rate of the gas, mg,;, can be expressed with respect to time, ¢, as shown in Eq. (1). 


om 
. — gas — 
Mgas = at. PgasAVgas (1) 


Equation (2) shows the mass flow rate balance between the inlet and exit of the control volume. The water vapor 
mass flow rate, my,, is treated as a mass source. 


PgasAVgas + Myy = (Pgas + OPgas (A + 0A) (Vgas + OVgas) (2) 


Equation (3) shows the gas momentum balance across the control volume. The expression below also approximates 
Vwy = Vgas, Stating that the velocity of the vaporized gas is approximately the speed of the flowing gas. 


(PgasAVgas)Vgas + MgasYwy + PA +(P + OP/2)0A 
= (Pgas + OPgas)(A + 0A)(Vgas + OVgas)(Vgas + OVgas) + (P + OP)(A + OA) (3) 


The gas conservation of energy expression becomes a bit more complicated. Equation (4) shows the general gas 
energy rate balance between the inlet and exit of the infinitesimal volume. A subscript of ‘in’ refers to an energy term 
entering the control volume, while a subscript of ‘out’ refers to an energy term as it exits the control volume. The 
control volume is otherwise adiabatic. In the Eq. (4) below, H is the enthalpy rate term, KE is the kinetic energy rate 
term, and Q-ony is the convective heat transfer term between the particles and the flowing gas. 


Basin + RE jase + Pye + RE pints + eee = Fog + RE rasmus an Beant + RE webu (4) 


The gas conservation of energy equation is expanded in Eq. (5) and also approximates Vy = Vgas. 


) + Did md? hi(T»i = Tyas) (#i) = 
(VgastOVgas)* 
a ee 


Dia . Pane 
PgasAVgas (aloes + ) + Myy (Gets + — 


2 
PgasAVgas (ome a OT jac) ae + Mwy Ge (eer 7 Ola) at (5) 


(VgastOVg “s)") 


The convective heat transfer term, Oa: in Eq. (4) is expressed in blue highlight in Eq. (5) as the sum of all the 
convective heat transferred from each individual particle size, i. In the above expression, 7 is the number of particle 
size bins and #; is the number of particles in the i” bin. 7, is the temperature of the particle, and / is the heat transfer 
coefficient for that particle. Also, Cgas and Cy in Eq. (5) are the heat capacity of the gas and water vapor, respectively. 

The speed of sound of a fluid, a, is expressed in Eq. (6), where k is the ratio of specific heats of a gas. 


2_ OP _ kp 


(6) 


OPgas - Pgas 


Equation (7) converts the conservation equations, which are normally temporally based, (i.e. 0/0f), into spatially 
based expressions (i.e. 0/Ox). The variable geometry along the tunnel axis, which is known, can be expressed as a 
differential equation, where v,eris the velocity of the reference frame. The reference frame for this model is the slowest 
moving particle velocity. The parameters in the tunnel are approximated to be steady state (i.e. not transient), and only 
vary in the axial direction. 


dA dAdt _ OA 1 4 
ax Ot Ax ~—«OE-« Vref (7) 


Combining and rearranging Eqs. (1), (2), (3), (6), and (7), while eliminating second order terms, provides 
expressions for the change in gas density, velocity and pressure across the infinitesimal volume and are given in 
Eqs. (8), (9) and (10) respectively. 
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OPgas _ PgasYgas(PgasAVgastMwyv) (5 as 1 wus) (8) 
ax (KPA—-Vgas(PgasAVgastMwy)) \AOX Vref Mgas 

Opis -kPAvgas (; aA 1 a) (0) 
ax (KPA—Vgas(PgasAVgastMwy)) \AOX = Vref Mgas 
oP = kPvgas(PgasAVgast+Mwv) Ges er we) (10) 
ax (KPA—-Vgas(PgasAVgastMwy)) \AOX Vref Mgas 


Combining, and rearranging Eqs. (1), (5) and (7), while eliminating second order terms, provides an expression 
for the change in gas temperature across the infinitesimal volume and is given in Eq. (11). 


. : 3 Ovgas 
OTgas _ PgasAVgasQconv—MwvPgasAVgasCwv(Tgas—Ty)—Vgas(PgasAVgasCgast+Mwv) J (11) 


Ox MgasV ref (PgasAVgasCgast+MwvCwr) 


Once again, the convective heat transfer term, Qcony, is the sum of all the individual heat transfers from all 
particles, as expressed and highlighted in blue in Eq. (5). 


C. Thermal Model Equations — Particle Conservation Equations 

The reader is encouraged to refer to the work by Bartkus et al.'? for the derivation of equations that describe 
changes in the particle mass, velocity, temperature and liquid content. These equations, shown below, were derived 
from the particle conservation of mass, momentum and energy equations. The appendix of that reference also provides 
all of the air, water, ice, and vapor thermo-physical properties as functions of temperature, pressure, and humidity that 
are used in this model. 

A mass balance is done for every particle size with the conservation of mass equation for particle size i shown in 
Eq. (12). In this equation, the subscript ‘p’ refers to the particle. Also, d refers to the particle diameter, while /», is the 
mass transfer coefficient. 


OMy j 1 
peels Jee 2 . = po 
ax mdi hm iPgas (yas Wy,i) Vref (12) 


Equation (13) uses the conservation of momentum to provide an expression for the change in particle velocity 
across an infintesimal distance for every particle size i. In Eq. (13), Cp refers to the coefficient of drag. 


2 
OVp 3 PgasCp,i(Vgas—Vp,i) 1 


= See es Pt (13) 


Ox 4 Ppidi Vref 
The energy rate balance provides an expression for the change in particle temperature for every particle size i as 


shown in Eq. (14). In this equation, 4 refers to the heat transfer coefficient, while Z is the latent heat of 
evaporation/condensation or sublimation/deposition. 


OMy j 
2h. -T.;)+L; 22 Dil 
OTyi _ 6(1d7 hi (Tyas-T pi) +24? hm, iPgasli(gas—@p,i)) = 6( ma? hi(Tyas Tpi)+Livres( ax )) 
ax Td PpiCpivres Td} PpiCpiVref 


(14) 
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For a particle undergoing phase change, the conservation of energy provides an expression, Eq. (15), for the change 
in particle melt fraction, 7, for every particle size i. In Eq. (15), Lmex refers to the latent heat of melting (1.e. latent heat 
of fusion). 


Om, ; 
2h. —T.; F meiees (7h) 
ONpi _ 6( Wd? Ni(Tgas-T pi) +24; NmiPgasli(@gas—Wp,i)) _ 6( nd? ni(Toas Tpi)tbivres Gx )) 


15 
Ox 1d} Pp ilmeltVref 1d} Pp ilmeltVref ( ) 

Since mass is conserved, the change in particle mass is the opposite change in the gas mass. Mass loss due to 
evaporation from a particle, for example, is a gain in gas mass. Water vapor is treated like a source term in the gas 
conservation equations, Eqs. (2) - (5), and is the sum of all the mass transfers for each particle as expressed in Eq. (16). 


OmMwy _ son [ee ) — Mwy 
ax i=1| ax eas ae (16) 


D. Instantaneous Condensation for Supersaturated Conditions 

At the time of writing, the authors were uncertain about what happens when conditions in the tunnel can potentially 
exceed saturation near the tunnel exit. For example, elevated humidity levels in the tunnel inlet (plenum) can lead to 
supersaturation at the tunnel exit, due to a rapid static temperature decrease from isentropic expansion in the tunnel 
contraction. Condensation occurs if conditions are supersaturated, but only at the rate as physically allowed through 
diffusion. With higher speed flows, it was not clear if the vapor potentially in excess of saturation would condense (or 
deposit as ice) quickly enough to not exceed saturation. 

For tests that may potentially exceed saturation, two simulations were run in an attempt to bound the problem, one 
where supersaturation was allowed, and one were the gas never exceeded saturation. For simulations that allow 
supersaturation, the equations as described in Sections II.B and II.C remain unchanged. For simulations that do not 
allow saturation to be exceeded, a subroutine was implemented that dealt with the mass that exceeded saturation which 
alters Eq. (12). 

In order to not exceed saturation, the algorithm implemented instantaneous condensation (IC). For purposes of this 
paper, instantaneous condensation is referred to the amount of vapor mass that must condense in order to not exceed 
saturation; condensation rates faster than can occur through diffusion alone. In the model, no new particles were 
generated from IC; all condensation occurred on existing cloud particles only. For modeling simplicity, the IC mass 
was equally distributed among all of the existing cloud particles. The mass rate form of this distribution is shown in 
Eq. (17). 


(17) 


Mpyi = Sn Hey 


In Eq. (17), yy is the vapor mass in excess of saturation, and is simply the m,,,, value in Eqs. (8) - (11). Eq. (18) 
shows the above equation in terms of a differential with respect to x. 


OMy j m 
Pl __ wv 
eS 18 
ax Vref Yeo (Hi) ( ) 
Equation (18) replaces Eq. (12) for calculating the changing particle mass. In addition, Equations (14) and (15) 
use this new mass increase as shown in the far right for each respective equation. 


E. Thermal Model Algorithm 

The thermal model was written in MATLAB version R2015b!’. The model solves the differential equations of 
mass, momentum and energy transfer. All conservation equations presented in the following model formulation 
section are solved using MATLAB’s built-in ODE45 solver. The system of differential equations is solved by first 
inputting the initial conditions. Time marching methods are employed by the MATLAB solver such that each iterative 
solution is solved with numerical relative and absolute tolerances of 10°. Physical accuracy of the solution is 
dependent on the accuracy of all the property values input into the model. For all simulations presented in this paper, 
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the energy transferred between the gas and particle(s) is balanced to the order of 10. The mass transfer between the 
gas and particle(s) initial and final states is balanced to the order of 10°). 


IfIl. Model Sample - Input And Output 


A. Experiment Descriptions 

Figure 2 depicts the PSL geometry used for the 2015 tests. The PSL icing tunnel had a 27:1 area contraction ratio 
from the plenum (Tunnel Inlet) to the test section (Tunnel Exit), with an axial distance of about 8.84 m. The spray 
nozzles and spray bar system are located at the tunnel inlet in the plenum, while the test section is a 0.91 m (3.0 ft) 
diameter free jet exit. PSL has the capability to spray either de-ionized water or city water (i.e. water that was not 
treated). City was was used for all tests conducted, with an approximate freezing temperature of 0 °C (32 °F). The 
injected cloud covered most of the cross-sectional area at any point in the tunnel, but radial variations in total water 
content (TWC) and particle size distribution (PSD) likely exist. During a given test point, the facility generally 
maintained target pressures, temperatures and plenum relative humidites to approximately +0.3 kPa (.05 psia), 
+ 0.5 °C (1 °F), and + 1%, respectively!. 


Inlet plane 


Free Jet 


Exit Plane 
(Test Section) 


Figure 2: PSL tunnel geometry used for 2015 tests with spray nozzles located at the tunnel inlet while the 
temperature, humidity and water content measuring instruments were located near the tunnel exit. 


Two experimental configurations are described for this paper and are shown in Figure 3. The first configuration 
was arranged to measure cloud content, gas temperature and humidity. Cloud content and phase were measured using 
the Science Engineering Associates (SEA) Multi-wire Probe!®”° which was located on the tunnel centerline near the 
exit of the free jet. Gas temperature and humidity measurements were made with a probe that was placed near the exit 
of the free jet, with the probe inlet offset 0.25 m (9.7 in) from the centerline. Struk et al.'> provides a description of 
the probe. In addition, all temperature, humidity and multi-wire data can be found in the same reference. Figure 3A 
shows the location of the multi-wire, and the combined temperature-humidity probe. A second configuration was 
arranged to measure particle sizes using a Cloud Droplet Probe (CDP) manufactured by Droplet Measurement 
Technologies, Inc. The probe was placed on the centerline axis near the exit of the tunnel and is pictured in Figure 3B. 
The CDP probe is capable of measuring particle sizes from 2 to 50 um. Van Zante et al.'° describes the configuration 
in greater detail and provides a summary of the cloud characterization. A description of the probe's operation is 
described by Van Zante and Rosine". 
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Multi-wire 


A Temperature and 
nelealteliavae)geler= 


Figure 3: Image of two configurations used at PSL where A) the multi-wire and the combined temperature-humidity 
probe are shown in the first configuration and B) the particle size measuring CDP probe is shown in the second 
configuration. Images are taken from the plenum camera, forward-looking aft. 


B. Model Input 

Experimental conditions are needed as inputs for the model. The main gas parameters required are the inlet gas 
temperature, pressure, and relative humidity, as well as the gas velocity at the tunnel exit. The main cloud inputs 
needed are the inlet median volumentric diameter (MVD), total water content, and cloud temperature. In addition, an 
approximate initial cloud velocity is required since this does not model the immediate exit of the pressurized spray 
nozzles. All sprays are initially liquid. Finally, the tunnel geometry, previously described, is the final model input. 

In order to more accurately simulate the cloud from the injection point just aft the spray nozzle bar system to the 
tunnel exit, an initial cloud PSD is required. The NASA PSL uses identical spray nozzles as in the NASA Icing 
Research Tunnel (IRT), and it is postulated that the measured PSD values at the NASA IRT test section are 
approximately equivalent to the initial PSD just aft of the nozzles at NASA PSL. PSD measurements were made in 
the IRT test section, but since the tunnel is a closed loop and near saturation, it is approximated that the PSD as it exits 
the IRT spray nozzles is the same PSD at the test section. PSDs of different cloud sizes were measured at the NASA 
IRT for various nozzle settings, namely the nozzle water pressure, Pyarer,noz, and nozzle atomizing air pressure Pair,noz”!. 
For this work, three initial PSDs were used, with median volumetric diameters, MVD, of 15, 20 and, 40 um 
respectively. Each distribution uses the PSD from the IRT and contains at least 24 bins. 


C. Model Sample Results 

Several parameter subscripts need to be defined, which will apply for the remainder of the document. Subscripts 
of ‘0’ and ‘s’ refer to a total or static value, respectively. A subscript of ‘i’ refers to a parameter at the tunnel inlet, 
while a subscript of ‘e’ refers to the tunnel exit. Conditions with a subscript of ‘off” refer to parameter values prior to 
the spray nozzles being activated, while a subscript of ‘on’ refers to values when the spray is turned on. A subscript 
of ‘exp’ refers to experimentally measured data, while ‘sim’ refers to model simulation results. 

Results from a simulation of Test Point # (TP#) 672 are provided as an example of the thermal interactions between 
the injected cloud and flowing gas. Figure 4 shows the predicted A) temperature, B) melt ratio, C) mass fraction that 
has evaporated, D) humidity along the tunnel axis, E) cloud residence time, and F) volumetric diameters. The inlet 
total temperature, 7); = 277.3 K (39.5 °F), inlet total relative humidity, RHo;= 10.8%, and inlet total pressure, 
Po; = 87.5 kPa (12.7 psi), were used as initial inputs into the model. A gas velocity at the tunnel exit of ve = 85 m/s 
was specified as a model parameter input. Liquid water was injected with an initial temperature of 280.4 K (45 °F) at 
a bulk total water content of TWCpun = 2.26 g/m>. The value of TWCoux was calculated as the ratio of water mass flow 
rate to volumetric gas flow rate at the test section and is calculated using Eq. (19) below. In this equation Myate, iS 
the water mass flow rate, and A is the cloud effective area which is assumed to be the same as the tunnel exit 
(i.e. 0.91m diameter). This 7WCp.% value assumes that all injected water reaches the exit of PSL (i.e. no evaporation 
or condensation) and is uniformly distributed across the exit area. TWCpuix is simply an input value whose value 
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changes as the cloud partially evaporates/sublimates before reaching the tunnel exit. The simulation used the 
MVD; = 40 um particle size distribution, a close representation of the approximate 42-um cloud that was initially 
injected into the tunnel by the spray nozzles. The PSD used in this simulation contained 44 bins. 


TWCOpuk = water (19) 
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Figure 4: Simulation of TP# 672 where A) temperature, B) melt ratio, C) mass fraction that has evaporated, D) humidiy, 
E) residence time, and F) volumetric diamters are shown with respect to the tunnel axial distance. 
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Figure 4A shows the temperatures of three differently sized particles in the cloud, as well as the gas static 
temperature, Ts sim, and static wet-bulb temperature, Twbs, sim, along the axis of the PSL icing tunnel. The three particle 
sizes correspond to the smallest (dia = 3 um), volumetric median (dia = 40 um), and largest particle (dia = 330 um) 
of the PSD used. The gas static temperature decreased from a combination of convective heat and mass transfer, as 
well as from isentropic expansion. With the cloud off, the gas static temperature at the tunnel exit calculates to be 
T:,e,off;sim = 273.7 K (33 °F), but with the cloud on, the gas temperature decreased to Ts,¢.onsim = 271.1 K (28.3 °F). The 
simulation predicted a 2.6 K (4.7 °F) decrease in gas temperature at the tunnel exit with the cloud activated. Also, all 
particles in the cloud decreased in temperature from the initial injection particle temperature. The largest particle was 
mixed phase which reached a temperature of 273.15 K (32 °F) at the test section. After completely glaciating, the 
40um-particle (and all other particles that froze out but are not shown in the figure) reached the static wet-bulb 
temperaure, Twb;. This can be seen in Figure 4A as the 40-um particle follows the calculated wet-bulb temperature 
curve (gray dashed line). Wet-bulb temperature is the temperature of a wet adiabatic surface undergoing evaporation 
(or sublimation of ice). 7wh, is the resulting temperature when convective heat transfer is balanced by evaporative 
cooling and is a function of 7;, RH, and P,*”. The smallest particle (3 zm) quickly froze out and evaporated/sublimated 
completely away, within traveling a distance of 0.037m in the tunnel. This is expected as the large surface area to 
volume ratio for small particles encourages heat and mass transfer and will have a fast response. 

Figure 4B shows the melt ratios, 7, for the three differently sized particles as well as the summed total melt ratio 
for the entire cloud. Melt ratio is defined as the ratio of liquid content to total water content. For a liquid water droplet, 
the melt ratio equals 1, while a solid ice particle, the melt ratio equals 0. In the simulation, the largest particle began 
to freeze, but remained mixed phase, and largely liquid, by the tunnel exit. This slow transition in phase is expected 
because the small surface area to volume ratio for large particles dampens heat and mass transfer and will respond 
more slowly. The 40-pm particle glaciated completely by 0.65 m down the tunnel. The smallest particle began to 
freeze at a distance of 0.001m, became completely glaciated by 0.004 m and as mentioned before, vaporized 
completely by 0.037m. Summing the liquid and total masses of all remaining particles, the final melt fraction at the 
tunnel exit totaled Nesim = 0.21. 

Figure 4C shows the mass fraction that vaporized for each of the 3 different particle sizes and the summed total 
amount that vaporized from the cloud. A value of | indicates complete evaporation or sublimation, while a value of 0 
means no mass change occurred. A negative value would indicate mass gain, through condensation or depostion. As 
expected, the largest particle experienced the least amount of mass loss due to its low surface area to volume ratio. 
Over half of the mass of the 40-um particle vaporized, while as mentioned before, the smallest particle vaporized and 
reached the value of 1 quickly, at which point it disappeared. Approximately 56% of the total cloud mass vaporized 
completely by the tunnel exit. 

Figure 4D shows the predicted mixed mass ratio, Wsim, (black dashed line) and static relative humidity, RH,sim, (red 
solid line) for the flowing gas along the tunnel axis. Mixed mass ratio is defined as the ratio of vapor mass to dry air 
mass. Mixed mass ratio increased along the tunnel axis from 0.64 to 1.83 g/kg, due to particle evaporation and 
sublimation. RHs,sim increased as well, from a combination of added vapor content to the flowing gas and decrease in 
static gas temperature. With the cloud off, the static relative humidity at the tunnel exit calculated to be 
RAs eoffsim = 13.3%, but increased to RAs,¢,on,sim = 47.6% with the cloud activated. 

The shaded region in Figure 4E represents the tunnel geometry and is shown for positional reference. The red 
curve in Fig. 4E represents the amount of time that has elapsed for cloud particles that had been sprayed out of the 
nozzles to reach a particular axial distance in the tunnel. It took cloud particles 1.25 seconds to travel from the spray 
nozzles to the exit of the tunnel. It is important to note that the cloud resided in the plenum for approximately 
1.0 seconds, or about 80% of the total time. The cloud particles resided in the plenum for so long because they were 
moving slowly at only a few meters per second in this section of the tunnel. The cloud particles accelerated through 
the contraction, and due to the elevated velocity, passed through the remaining tunnel section quickly. The largest 
changes in total gas temperature, melt fraction (freezing), mixed mass ratio (absolute humidity), and particle 
evaporation occurred in the plenum, since the cloud particles resided here for a substantial fraction of the time. 

In an attempt to characterize the cloud particle size distribution, Figure 4F shows three volumetric diameters along 
the axis of the tunnel, Dv/0, Dv50, and D90. Volumetric diameter is a representative diameter where the specified 
percentage of the total volume of the spray is composed of particles with diameters smaller than or equal to the stated 
value. For example, Dv10 is the diameter where 10% of the cumulative volume of the spray is made up of particles 
smaller or equal to the specified diameter. Dv50, which is by definition the MVD, signifies the diameter where half of 
the cumulative volume is composed of particles smaller than this diameter while the other half of the total volume is 
greater than this diameter. Finally, Dv90 represents the diameter where 90% of the cumulative volume is smaller than 
the diameter calculated for Dv90. The Dv/0 and Dv90 values provide a width of distribution, or a band around the 
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median Dv50 value. All three Dv curves increased from their initial value, despite cloud mass loss due to 
evaporation/sublimation. In this example, Dv50 increased from 401m to 821m from the tunnel inlet to the exit. This 
increase in Dv is caused by smaller particles completely evaporating, leaving only the larger particles. Smaller particles 
have a greater surface area to volume ratio and evaporated more quickly than larger particles. Even though these 
remaining large particles were smaller than their initial size, the mathematical calculation of the three Dv values 
increased for this time period. It should be noted that provided enough time and the proper conditions to evaporate 
even the largest particles, the Dv curves would eventually all decrease to 0 um. 


IV. Model — Experiment Comparisons and Discussions 


The focus of the May 2015 testing was to examine the effect various PSL facility parameters have on the melt 
ratio of the cloud at the test section. Three primary parameters were examined: bulk TWC, plenum humidity, and 
spray bar water temperature. Within these sets of results, target particle size was also varied between a smaller, 
nominally 15-um MVD);, and larger value of 40-um MVD;. Model results for each of these test variations are discussed 
in separate sections below. Model results are also compared with experimentally measured particle size data in a 
separate section. 


A. Total Water Content Sweeps and Comparisons 

For the rest of the document, any references to temperature in refer to the gas temperature, and the subscript of 
‘gas’ is omitted for ease of reading. Table 1 below shows the results of two sets of tests where bulk total water content, 
TWCouk was varied. Facility target conditions are presented at the top half of the table for each test set, which includes 
the inlet total pressure, Pj, static pressure at the tunnel exit, Ps.o9; gas velocity at the test section, V., inlet total 
temperature, 7), and static temperature at the tunnel exit, Ts,¢,o4: The value of Po, is approximated to be constant for 
the entire length of the tunnel and is independent of spray on or spray off conditions. Also, Ps,co4 is approximated to 
be independent of spray on or spray off condition at the tunnel exit, and can also be written as Ps. The relative 
humidity in the icing tunnel plenum was specified during testing, which is a close approximation for the total relative 
humidity, RH, in the slow moving plenum. The corresponding mixed mass ratio in the plenum, i, is also provided. 
From the temperature, humidity and pressure values, the total wet-bulb temperature, Twbo,io#, and the static wet-bulb 
temperature, Twhs,e,o4 are provided. Wet-bulb temperatures were calculated using equations as defined in Veres et 
al.”?. The above parameters are facility target conditions in the absence of the icing cloud. The initial cloud 
characteristics (immediately aft the spray nozzles) are given in the table, which include the bulk total water content, 
TWCoux initial water temperature, Tiwaer,, and initial cloud median volumetric diameter, MVD;. The MVD; was 
calculated using specified nozzle conditions, namely the nozzle water pressure, Piyater,noz, and nozzle atomizing air 
pressure Pair noc’). 

The bottom half of the table compares model results with experimentally measured values. The focus is on the 
change in experimental conditions when a cloud is present. The change in mixed mass ratio between cloud on and off 
was experimentally measured, AWex,. The corresponding simulation value, AWsim, is given, as well as how closely the 
simulated value matches the experimentally measured value, Awaz. The value of Awaz is calculated as the 
following: 


AW sim- AWex 
AW aif f = esta eee x 100% (20) 


AWexp 


The change in gas total temperature at the test section between cloud on and off, ATo,¢,exp, is provided, along with 
the corresponding simulated value, AT .,sim. Similarly, the measure of how closely the simulated value matches the 
experimental value is given as AT ,¢ “aig, and is calculated as follows: 


ATo,e,sim— ATo,e, 
AT ¢, dif f = eaieeese cexP x 100% (21) 


In Egs. (20) and (21) a negative value means the simulation under-predicted the experimentally measured change, 
while a positive value means the model over-predicted the change. 

The melt ratio is defined as the mass ratio of liquid water content to total water content. The experimentally 
measured melt ratio, 77,cxp, and simulation value, 7,sim, at the tunnel exit are provided. The 7.x) value was calculated 
by taking the maximum value between the two LWC wires of the multi-wire (the 2.1 mm and 0.5 mm wires) and 
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dividing by the measured value from the TWC wire (the half-pipe)!°. No correction factors were applied, as only the 
ratio of the raw measurements from the multi-wire are reported. Finally, the change in total wet-bulb temperature at 
the test section between cloud on and off, ATwhoeexp, iS given, as well as the corresponding simulation value, 
ATwhoe.sim- 


Table 1: Total Water Content Data Set 1 Sweeps 
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Between the two data sets, experiments were conducted at a higher ambient pressure in the first data set, TWC Sweep 
1, which also meant higher wet-bulb temperatures. The initial cloud size for each test was MVD; = 15um. In both 
sweeps, the experimentally measured vapor content increased at the test section as TWCp.% increased. Simulations 
matched this humidity trend within approximately 30% of experimentally measured change for all 7 tests. As TWCouix 
increased, the drop in gas temperature generally increased. The model similarly matched this temperature trend to 
within approximately 30% of the measured change for all 7 tests. Comparing similar conditions from the two data 
sets, namely TP# 666, 667, and 665 vs TP# 678, 679, and 680 respectively, it can be seen that the increase in humidity 
was greater for the lower ambient pressure tests. Similarly, the decrease in temperature was greater for the lower 
pressure data set. These trends were expected and predicted by simulation as there is greater potential for evaporation 
at lower pressures. 
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For TWC Sweep 1, all simulations predicted a fully glaciated cloud. For the two higher TWCp,% tests of 
TWC Sweep | (TP# 665 and 664), the 7, cx) was measured as 0.10. Such a low value suggests a mostly or fully glaciated 
cloud since a fully glaciated cloud produces a signal (~5% of IWC) on the LWC elements of the multi-wire resulting 
in a non-zero 7. For the two lower TWC},x tests, in particular TP# 666, higher 77¢,cx, values were measured, suggesting 
a mixed phase cloud. All things being equal, it is not expected to see a more liquid cloud as the TWC},1% decreases, 
and in fact, the 7c.) measurement is expected to be more glaciated as TWCbi decreases. It is hypothesized that the 
resulting particle size distribution, for TP# 666 in particular, may have had a heavier right tail particle size distribution 
than what was used for simulation. While according the IRT particle size distribution calculator?! the nozzle settings 
produce a 15-pm cloud, the atomizing air pressure of Pairnoz = 103 kPa (15 psig) was lower for TP# 666 than the rest 
in the data set. As a result, the choked flow in the nozzle produced a weaker resulting shock , or expansion fan waves”?, 
as the gas expanded out of the nozzle compared to the greater pressurized nozzle tests. This weaker shock may have 
resulted in a greater number of larger particles, which would not have frozen out as quickly as smaller particles. In 
addition, there is a tendency for larger particles to concentrate towards the center of the tunnel, due to a combination 
of particle inertia and tunnel contraction. The multi-wire was measuring the melt ratio values at precisely the tunnel 
center and may also have contributed to this higher 7%,cx, measurement for TP# 666. Particle size measurements were 
not made for these tests or any of the parameter sweeps, so as to help verify this hypothesis. 

Similarly, all three simulations predicted full glaciation for TWC Sweep 2. While TP# 680 suggests a mostly or 
fully glaciated cloud, the two lowest TWCiui tests, TP# 678 and 679, measured mixed phase clouds. Again, the 
atomizing air pressure, in particular TP# 678 was low, Pairnoz = 103 kPa (15 psig), which when expanded out of the 
nozzle created a weaker shock that generated some large particles that may not have fully glaciated. This may account 
for the mixed phase experimental value, even though it is not expected for this low TWCp,, case. One thing to note is 
that the experimentally measured melt ratio of TP# 678 (Nery = 0.27) is lower than its similar counterpart TP# 666 
(Neexp = 0.79). Since TP# 678 was conducted in a lower ambient pressure, more evaporative cooling occurred, which 
reduced the melt ratio slightly more than the higher pressure companion test. 

A trend is recognizable in wet-bulb temperature, Twb, as well. As TWCp,i increased, an increase in ATWHo,¢,exp WaS 
measured, which was also predicted by simulation. The increase was not large as there was a near balance in the 
calculated Twb using the pre-spray humidity and temperature and the calculated 7wb value during spray-on where 
there was an increase in vapor and decrease in gas temperature from evaporating particles. This increase in ATwbo is 
attributed to the energy added to the system from an initial water temperature that was greater than the wet-bulb 
temperature. An assertion made for calculating wet-bulb temperature is that the evaporating water temperature is 
supplied at the wet-bulb temperature. For an unsaturated, adiabatic system where the water temperature is greater than 
Twb, the calculated wet-bulb temperature will increase and the final Twb value will be dependent on the amount of 
water mass (or more specifically, the amount of energy associated with the water mass) in the system. Hence, higher 
TWCobux experiments in Table | provided more energy to the gas since Twarer > Twp for every test, which resulted in 
greater ATwho,- values compared to lower TWCy.% tests for a given TWC sweep. 

In addition to the change in wet-bulb temperature, information on the state of the cloud can be gleaned from the 
absolute 7wb. For example, the total and static wet-bulb temperatures for TWC Sweep 1, Twhoioy = -2.0 °C (28 °F) 
and Twby,eof = -4.2 °C (24 °F) respectively, were both below freezing. This means that once injected into the tunnel 
plenum, the evaporating particles experienced sub-freezing conditions, even with an approximate 0.5 °C (~1 °F) 
increase in total wet-bulb temperature. With a freezing point of 0 °C (32 °F) for city water, it is expected to see a fair 
amount of glaciation, especially for clouds composed of smaller particles such as in these two data sets. 
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Table 2: Total Water Content Data Set 2 Sweeps 
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Two series of TWCp.% sweeps are presented in Table 2, where the target total temperature of TWC Sweep 3 is about 
2.4°C (4°F) warmer than TWC Sweep 4. The initial cloud MVD; for each test was 40 um. As was the case before, the 
experimentally measured vapor content increased as TWCp, increased. Simulations matched this humidity trend 
within approximately 20% of experimentally measured change for all 8 tests. Again, as TWCp.,% increased, the drop 
in gas temperature generally increased. The model similarly matched this temperature trend to within approximately 
40% of the measured change for all 8 tests. Comparing similar conditions from the two data sets in Table 2, it can be 
seen that the increase in humidity was slightly greater for the warmer test set. Similarly, the decrease in temperature 
was slightly greater for the warmer data set. These trends were expected and predicted by simulation as there is greater 
potential for evaporation at warmer temperatures. The differences were small perhaps due to the fact that the 2.4 °C 
(4 °F) between the two data sets in Table 2 is not a large difference. 

Simulations predicted a melt ratio of approximately 7,sim = 0.20 for all 4 tests runs in TWC Sweep 3, regardless of 
TWCouk. This value matches very closely the melt ratio measured for the for the two higher TWC tests, TP# 672 and 
TP# 673, where 7,exy = 0.23 and H,exp = 0.27 respectively. As was the case before, a lower TWCopuix should not result 
in a more liquid cloud, but this was the case for TP# 670 and TP# 671 where 7, exp = 0.69 and 7.,exp = 0.66 respectively. 
Again, it is speculated that the lower nozzle air pressure settings for the lower TWCoux tests produced a heavier right 
tail PSD than what was used for simulation. The same trend appears again for TWC Sweep 4. The model predicts quite 
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well the melt ratio for the higher TWCyu. tests, but not for the lower TWCp.% tests, again perhaps due to the lower 
nozzle air pressure setting. It can be seen that the melt ratios for the warmer test are higher, which is expected and 
predicted by the model. 

Despite sub-freezing Twhoiog¢ and Twhseog values for both data sets, the cloud did not glaciate completely, 
experimentally or in simulation. The 40-um MVD; clouds contained some large particles that did not fully freeze before 
reaching the tunnel exit. Also, there again exists a positive correlation between the increase in ATwhoeexy and 
increasing TWCp.%, which is captured in simulation. As mentioned before, this trend exists because of the additional 
energy the gas receives from water temperatures in excess of the wet-bulb temperature. 


B. Relative Humidity Sweeps and Comparisons 

Three tests were run that varied RHo,; and are listed in the RH Sweep / data set in Table 3. The tests were run 
simulating a low altitude with TWCs.x% = 1.0 g/m?. As RHp,; increased, the magnitude of A@e,ex» decreased as did the 
magnitude of the gas temperature drop. This trend is expected as higher RH; values reduce the potential for 
evaporation, and these trends were captured in the simulation to within 30% of the measured value changes. 


Table 3: Relative Humidity Data Set Sweeps 
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For TP# 667, where RHp; = 10%, a value of 7%,» = 0.16 was measured, suggesting a mixed phase, but perhaps a 
nearly glaciated cloud. The simulation predicted a fully glaciated cloud. With below freezing wet-bulb temperatures, 
Twhoiof = -2.2 °C (28 °F) and Twhs, coy = -4.4 °C (24 °F), a mostly glaciated cloud is expected when the cloud MVD; 
is small (15 um). When RHo,; was increased to 35% for TP# 668, a value of 7,cx) = 0.81 was measured, suggesting a 
nearly all liquid cloud. The simulation predicted nearly the same value with 77. sim = 0.82. Since the cloud spent most 
of the time in the slow moving plenum experiencing above freezing conditions with a total wet-bulb temperature of 
approximately Twhoiog= 0.6 °C (33 °F), and then for a brief moment static wet-bulb temperatures as low as Twhs cof 
= -1.7 °C (29 °F) during the acceleration period, it is not expected that much of the cloud would have glaciated. As 
RHpo,; was increased to 50% for TP# 669, 7,ex» was measured to be 0.82, not very different from the previous test, and 
again suggesting a mostly liquid or perhaps fully liquid cloud. The simulation predicted a fully liquid cloud. Again, a 
nearly all liquid cloud is expected since Twhoiog = 1.7 °C (35 °F), for which the cloud experiences for much of the 
duration in the icing tunnel, and for a brief moment experiences a barely freezing static wet-bulb temperature of about 
Twbs,e,of = -0.6 °C (31 °F) as it exited the tunnel. It is possible that the cloud remained fully liquid as it traveled down 
the length of the wind tunnel for TP# 669 as wet-bulb temperatures tend to increase slightly when the spray cloud is 
activated. 

A second data set was run where RHo,; was varied, this time at a higher altitude and higher total water content 
(TWCux = 4.9 g/m*). Simulations of this data set was run in two different ways. The first method allowed the flowing 
gas to become supersaturated and that data set is labeled RH Sweep 2 (RHsim > 100%). The second method limits the 
flowing gas to never exceed saturation (RHsim <= 100%) and is labeled as RH Sweep 2 Saturated in Table 3. In the 
latter simulations, instantaneous condensation is implemented into the algorithm if conditions approach 
supersaturation so as to not exceed saturation. In the former simulations, supersaturation is considered because it may 
be possible that as the gas reached near saturation levels in the plenum, and then accelerated in the contraction, the 
gas static temperature rapidly decreased via isentropic expansion with a mostly unchanged absolute humidity, 
resulting in possibly supersaturated conditions. The two methods are used to simulate this data set, in an attempt to 
bound the conditions. 

In the RH Sweep 2 data set, all three simulations exceeded saturation levels at the tunnel exit, with all three values 
around RH%s,¢0n,sim = 120%. The simulations in this data set over-predicted the absolute humidity increases and gas 
temperature decreases by about 30% of the measured values. When simulations were run not allowing supersaturation 
(RH Sweep 2 Saturated data set), the predicted absolute humidity increases and gas temperature decreases were 
predicted to within 10% of the experimentally measured values. This improvement in prediction suggests that near 
instantaneous condensation does exist. 

The predicted melt ratio matched the measured values quite closely for both data sets. This means that the melt 
ratio did not change much between the supersaturated and saturated simulations. Again one is able to make general 
predictions of the melt ratio looking at just the total and static wet-bulb temperatures. For TP# 681 and 682, both wet- 
bulb temperatures were sub-freezing, so it is expected that a small particle cloud would largely glaciate. For TP# 683, 
where the cloud experienced near freezing conditions in the plenum at Twhoo7= 0 °C (32 °F), and then very briefly 
experienced at its coldest Twhs,coy = -1.7 °C (29 °F) during the acceleration portion, it is expected that the particles 
would be mixed phase. It is also worth pointing out that the cloud went from nearly fully glaciated for TP# 682 where 
RHoi = 25%, to nearly fully liquid for TP# 683 where RH, = 50%. It is expected that the window to generate a mixed 
phase cloud using relative humidity as the parameter to control 7 is smaller for small MVD; clouds and larger for large 
MVD; clouds. 


C. Water Temperature Sweeps and Comparisons 

Two series of Tyaier,; Sweeps are presented in Table 4, where Twater Sweep I tests were conducted with an initial 
cloud size of MVD;= 15 pm, and Twater Sweep 2 was conducted with a cloud of MVD; = 40 um. Initial water 
temperatures were increased to as much as 82 °C (180 °F). The change in humidity was relatively insensitive to the 
increasing Tyaer,i for both sweeps. The model also predicted this nearly constant absolute humidity increase to within 
approximately 30% of the measured values. In addition, the drop in gas temperature was also insensitive to the 
increasing Twater,;, With the model predicting this insensitivity to within approximately 40% of the measured values. 
The model actually predicted slightly smaller decreases in gas temperature as nozzle temperature increased, due to the 
fact that the initially warm cloud transferred its heat to the gas, raising the gas sensible heat. These elevated water 
temperatures (temperatures greater than Tw) also account for the increasing ATwhoesim values as Tyare,i increased. 
The small 7TWCp.% values is the reason for the small A7Twbo.¢sim increase and the negligible differences in 
experimentally measured values of ATwho,e,exp aS Twater,i increased. 
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Table 4: Water Temperature Data Set Sweeps 
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FeO : 


3 
589 


8 
1 
1 


5 
6 
0 
ae 
0 
15 


<— 
| 0.69 | 

AT 0,0, aig (Yo) 
0.83 
Pesin() | 0.00 


The model predicted a fully glaciated cloud for all initial water temperatures in the Twater Sweep J data series, 
even for the initially high water temperatures for TP# 690 and 691. Rather high melt ratios, however, were 
experimentally measured for these two tests. The model predicted mixed phase conditions for TP# 689 and 692 in the 
Twater Sweep 2 data set, but predicted more glaciated values than what was measured. Further testing and analysis 
are required to understand this melt ratio discrepancy with respect to warmer Tater; values. 


D. Cloud Characterization — Particle Size Comparisons 

A total of 11 experiments where particle size data was collected is used to compare with model results. All 11 tests 
produced clouds with a presumed initial size of MVD; = 15 um. Only these smallest particle clouds were capable of 
being measured by the CDP probe alone as the PSD for these clouds are largely contained within the 2-50 um range 
of the probe. The other main experimental conditions varied for these 11 experiments with the minimum and maximum 
value for each parameter listed in Table 5. 


18 


American Institute of Aeronautics and Astronautics 


Table 5: Cloud Characterization Test Parameter Ranges 


All 11 experiments were run at conditions such that the model predicted potential supersaturation. Therefore, to 
bound the possible cloud size at the tunnel exit, MVD., two simulations were run for each experimental test, one where 
the tunnel gas was allowed to supersaturate, the second where the gas was limited to never exceed saturation. Figure 
5 isa 1:1 comparison graph of the predicted MVDe sim versus the measured MVD¢ exp at the tunnel exit. The red triangles 
in Figure 5 represent simulations that allowed supersaturation, while the blue circles represent simulations that limited 
the gas to never exceed saturation. The closer a simulation result matches the experimentally measured MVDe exp, the 
closer the marker will be to the 1:1 line. 
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Figure 5: A 1:1 comparison graph of predicted MVD, values with measured MVD. values. 


In general, the change in particle size was not large, according to experimentally measured values and simulation 
predictions. With an approximate initial experimental value of MVD; = 15 um, the measured final values increased to 
somewhere between MVDe,exy =16 and 20 um. The predicted MVDesim values for the supersaturated simulations also 
fall into this 16 to 20 um range, and generally fall on or above the 1:1 line. The predicted MVD¢,sim values for the 
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saturated simulations range from 13 to 17 um, and generally fall on or below the 1:1 line. On average, the MVDe.sim 
values for a supersaturated simulation was about 3 pm larger than saturated simulations. 

For the simulations where the gas became supersaturated, the cloud evaporated to the point where the smallest 
particles completely vaporized, leaving just the largest particles behind. The larger particles also become smaller from 
evaporative mass losses, but the overall effect mathematically resulted in a larger MVDesim values. In addition, 
particles that glaciated expanded in volume, which also played a role in increasing the cloud MVDesim value at the 
tunnel exit. 

The conditions for all 11 experiments, even prior to turning the spray nozzles on, resulted in a potentially 
supersaturated gas at the tunnel exit due to an isentropic reduction in static temperature. As a result, the cloud would 
have netted a gain in water mass if the gas was forced to instantaneously condense (IC) so as to never exceed saturation 
as the cloud moved through the contraction. The IC algorithm distributed the IC mass equally among all of the existing 
particles. Unlike when evaporation and sublimation increased MVDesim, condensation and deposition decreased 
MVDesim. While all particles increased in volume, the larger share of volume was gained by the smaller particles, 
mathematically shifting MVDe,sim values below the initial value. 

Overall, the model does a fair job in predicting final MVD, when the initial MVD; is small. The supersaturated and 
saturated simulations bound the experimentally measured values to within a few microns. Clouds composed of larger 
particles were run and measured during the most recent fundamental physics of ice-crystal icing tests at PSL in March 
2016. Further simulations are anticipated to be run with these larger MVD; clouds and compared to the experimental 
data. 


V. Conclusion 


The ultimate goal of the model is to better understand the complex interactions between the test parameters and 
have greater confidence in the conditions at the test section of the PSL tunnel. The model was developed in an attempt 
to explain the observed changes in test conditions for PSL icing tests by coupling the conservation of mass and energy 
equations for both the cloud particles and flowing gas mass. When compared to experimental data, the model predicted 
to within 30% of the measured gas temperature and humidity changes for a variety of parameters. This included 
variations of total water content, plenum relative humidity, and spray bar water temperatures. A subroutine was 
implemented that allowed for the option to control humidity so as to not exceed saturation at any point. Model results 
showed better agreement when supersaturation was not allowed. The model did reasonably well in predicting melt 
content at the test section, especially for clouds with larger particle sizes. There was some disagreement in melt 
fraction when initial water temperature was elevated to 43 °C (110 °F) and 82 °C (180 °F). Further testing and analysis 
is required to determine this discrepancy. In addition, the model predicted particle size at the tunnel exit with good 
agreement, however, the comparison was limited to clouds of an initial particle size distribution of MVD; = 15 um. 

One of the key findings from this work is that there was a nearly constant but slight increase in wet-bulb 
temperature when the spray cloud was activated for every test and simulation. This increase is wet-bulb temperature 
is attributed to the initial water temperature being warmer than than the wet-bulb temperature. Another finding from 
this work was that the total wet-bulb temperature in the plenum was a large factor in determining cloud phase. Finally, 
the model aided the fundamental physics of ice-crystal ice accretion tests in March 2016 in guiding the development 
of the test matrix and helping to identify the mixed-phase parameter space. 
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